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Temporal relaxation of density fluctuations in supercooled liquids near the glass transition occurs 
in multiple steps. The short-time /3-relaxation is generally attributed to spatially local processes 
involving the rattling motion of a particle in the transient cage formed by its neighbors. Using 
molecular dynamics simulations for three model glass-forming liquids, we show that the /3-relaxation 
is actually cooperative in nature. Using finite-size scaling analysis, we extract a growing length-scale 
associated with /3-relaxation from the observed dependence of the /3-relaxation time on the system 
size. Remarkably, the temperature dependence of this length scale is found to be the same as that 
of the length scale that describes the spatial heterogeneity of local dynamics in the long-time a- 
relaxation regime. These results show that the conventional interpretation of /3-relaxation as a local 
process is too simplified and provide a clear connection between short-time dynamics and long-time 
structural relaxation in glass-forming liquids. 


I. INTRODUCTION 


Temporal relaxation of density fluctuations in su¬ 
percooled liquids near the glass transition occurs in 
multiple steps. At short times, the temporal auto¬ 
correlation of density fluctuations and related corre¬ 
lation functions approach a plateau after a fast ini¬ 
tial decay. This part of the relaxation is known as 
the /3 relaxation. The subsequent long-time decay of 
the correlation functions from the plateau to zero is 
known as the a relaxation. The origin of this non¬ 
exponential, multi-step relaxation and possible con¬ 
nections between the short-time (3 relaxation and the 
long-time a relaxation are among the fundamental 
issues in the study of glassy dynamics. 



Many approaches to understanding slow relaxation 
in glass forming liquids invoke the notion of a grow¬ 
ing length scale that governs the increase of relax¬ 
ation time scales [2-§. Experimental and theoretical 
studies of dynamic heterogeneity in glass-forming liq¬ 
uids (M3 have lead to a detailed analysis of length 
scales that are associated with spatial correlations of 
the mobility of particles ITdMTsi l . These correlations 
have been studied through a four-point density corre¬ 
lation function <74 (r, t ) [191-121] . its Fourier transform - 
the four-point structure factor £4 (< 7 , t) - and the asso¬ 
ciated dynamic susceptibility = lim^o (S 4 (< 7 , t) 

[ 8 , H, [ 3 , HU, [23}. Analytic predictions for the be¬ 
haviour of dynamic length scales and susceptibili¬ 
ties in both short-time (/3) and long-time (a) re¬ 
laxation regimes have been obtained from inhomo¬ 
geneous mode coupling theory (IMCT) 0 0 EH- 
These predictions include an initial power-law growth 
of the cooperativity length scale in time in the /? relax¬ 
ation regime, followed by a saturation at time scales 
comparable to the a relaxation time. However, details 
of how the dynamics crosses over from its short-time 
behaviour to that at long times and the correspond¬ 


FIG. 1. Four-point susceptibility X 4 (t) for different sys¬ 
tem sizes, plotted versus time t. Note the appearance of 
a peak at short times. The height of the peak increases 
with increasing system size. The filled symbols represent 
data from molecular dynamics simulations and open sym¬ 
bols show data from Monte Carlo Simulations. The peak 
of X 4 at short time scales disappears in Monte Carlo sim¬ 
ulations. Inset: The height of the short-time peak of \4 is 
plotted as a function of the length L of the sample. The 
height is found to increase linearly with L. 


ing crossover in the dynamic length scale(s) are still 
not well-understood. There are many ideas and ob¬ 
servations that attempt to relate dynamical features 
observed at short time scales to long-time structural 
relaxation [24l43l( | . Thus, a proper understanding of 
relaxation processes in both /3 and a regimes and their 
underlying relation (if any) is extremely important in 
the overall understanding of glassy dynamics and its 
rich phenomenology. 

In an earlier study [U, the method of finite-size 
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scaling (FSS) [13, used extensively for obtaining ac¬ 
curate numerical results for critical properties near 
conventional phase transitions, was used to obtain the 
length scale £(T) associated with dynamical hetero¬ 
geneity, by analyzing the size dependence of the dy¬ 
namic susceptibility Xi{t) and the associated Binder 
cumulant. It was subsequently shown that the hetero¬ 
geneity length obtained in this way is in good agree¬ 
ment with that obtained by analyzing the four-point 
structure factor Si(q,t) [35j,|36j. However, very large 
system sizes needed to be employed in the analy¬ 
sis of Si(q,t), and further difficulties arise from the 
ensemble dependence of the associated susceptibility 
Xi(t) dH. Indeed, results HU claiming to confirm 
IMCT predictions for Xi{t) and the associated length 
scale £(T) are affected by finite-size effects. Thus, 
notwithstanding caveats about its use [13, Eo|, FSS 
offers an attractive approach to study length scales 
relevant to glassy dynamics. We employ this method 
in the present work. 

It was found, however, in 34] that the a relaxation 
time r a that describes the long-time decay of a two- 
point density correlation function does not exhibit the 
expected dynamical scaling, in that (a) the time scale 
decreases instead of increasing with system size, and 
(b) a plot of r a scaled by its asymptotic value for large 
system sizes vs. similarly scaled values of X4 does 
not display a scaling collapse. This was interpreted 
as indicating a mixing of activated and non-activated 
mechanisms of structural relaxation (see also [4lU42| ). 
It has been suggested [13 that behaviour described by 
IMCT, without the influence of activated dynamics, 
should be expected instead in the short-time dynam¬ 
ics, in a temperature window above the mode cou¬ 
pling transition temperature. 



FIG. 2. Mean square displacement (MSD) and the deriva¬ 
tive of the logarithm of MSD with respect to the logarithm 
of time f, shown as a function of t. The minimum of the 
derivative defines rp. (see text for details). 


In this article, we describe the results of FSS stud- 



FIG. 3. Top panel: rp for the 3dKALJ model, shown 
as a function of system size for different temperatures in 
the interval T G [0.45,0.90]. rp increases with increasing 
system size, saturating at a value that increases with de¬ 
creasing temperature. Bottom panel: Scaling collapse of 
rp. For each temperature, the linear size L is scaled by a 
length £(T) such that data for all temperatures collapse 
to a master curve. Inset: Comparison of the correlation 
length extracted from the system-size dependence of rp 
with that obtained from the a regime [HJ as a function of 
temperature. The two length scales agree with each other 
to a good accuracy. 


ies of a short time scale (the /? relaxation time scale 
Tp) for three generic model glass-formers. The simu¬ 
lation details are provided in Sec. II. We compare the 
estimated length scales with those obtained from FSS 
analysis of Xai the peak value of X 4 (t), which quan¬ 
tifies dynamic heterogeneity at the a relaxation time 
scale, and find that the two length scales agree well 
for all the systems we study. We also find a power-law 
dependence of the /? time scale on the length scale, in 
qualitative agreement with the prediction of IMCT. 
These results are described in Sec III. We also per¬ 
form an analysis of the dependence of the behaviour 
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FIG. 4. Top left panel: System-size dependence of > the peak value of for different temperature for the 3dIPL 

model system. Top right panel: The data collapse of Xa (-W T) to extract the dynamic heterogeneity length scale. 
Xi {N -A- oo, T ) is the asymptotic value of x± i n the limit of infinite system size. The corresponding length scales are 
shown in the inset of the bottom right panel of this figure (red circles). Bottom left panel: System-size dependence of rp 
for all the temperature studied, including high temperatures. Bottom right panel: Data collapse of rp (N, T ) to obtain 
the length scale associated with the cooperativity in /3-relaxation. Inset: Comparison of the correlation length extracted 
from the system-size dependence of rp with that obtained from the system-size dependence of Xa > as a function of 
temperature. The two length scales agree with each other to a good accuracy. 


in the /3 regime on the microscopic dynamics, as de¬ 
scribed in Sec. IV. Sec. V contains a discussion of 
our results and a summary. 


II. METHOD AND SIMULATION DETAILS 


tern characterized by a repulsive inverse power-law 
potential (3dR10) |46j whose range of interaction is 
smaller compared to that in the first two models. 
For the first two model systems the interaction range 
covers the second neighbouring shell whereas for the 
third case it is truncated within the first shell. 


We study three model liquids in this work. They 
are: (i) the three-dimensional Kob-Andersen binary 
Lennard-Jones mixture (referred here as 3dKALJ) 
[UJ, (ii) a three-dimensional system characterized by 
a repulsive inverse power-law potential (3dIPL) con¬ 
sidered in Ref. [45| and (iii) a three-dimensional sys- 


We have performed a series of simulations, for sys¬ 
tems of sizes ranging from N = 150 to N = 28160 par¬ 
ticles at number density p = 1.20 of the glass-forming 
Kob-Andersen binary liquid mixture 4J|. Details of 
the model parameters and the reduced units used 
for length, time and temperature are the same as in 
13411. We have done simulations for 7 different tern- 
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FIG. 5. Top left panel: System-size dependence of Xa , the peak value of X 4 (t), f° r different temperature for the 3dR10 
model system. Top right panel: The data collapse of Xa (-W T) to extract the dynamic heterogeneity length scale. 
Xa {N — > oo, T) is the asymptotic value of Xa i n the limit of infinite system size. The corresponding length scales are 
shown in the inset of the bottom right panel of this figure (red circles). Bottom left panel: System-size dependence of rp 
for all the temperature studied, including high temperatures. Bottom right panel: Data collapse of rp (N, T ) to obtain 
the length scale associated with the cooperativity in /3-relaxation. Inset: Comparison of the correlation length extracted 
from the system-size dependence of rp with that obtained from the system-size dependence of Xa j as a function of 
temperature. The two length scales agree with each other to a good accuracy. 


peratures in the range T £ [0.900,0.450]. For the 
3dIPL model the system sizes studies are in the range 
N £ [150,10000] at number density p = 1.20 for 9 
different temperatures in the range T £ [0.450,1.000]. 
For the 3dR10 model system sizes studied are in range 
N £ [64,4000] at number density p = 0.81 for 9 
temperatures in the range T £ [0.50,1.000]. All the 
simulations are done in the canonical ensemble us¬ 
ing a modified leap-frog integration scheme with the 
Berendsen thermostat. We have also performed simu¬ 
lation with another constant temperature simulation 
algorithm due to Brown and Clark j47|. The results 
do not depend on the exact algorithm used for inte¬ 
grating the equations of motion. The equilibration 
runs for all these temperatures are close to 100 r 


(where r is the a-relaxation time estimated from the 
decay of the two-point density correlation function) 
and we have averaged the data over 32 independent 
runs of length 100 r. For the results discussed in 
Sec. IV, we performed Brownian Dynamics simula¬ 
tions using the predictor-corrector algorithm given in 

la¬ 
in (34], the a relaxation time was identified by the 
decay of a two-point density correlation function, as 
well as by considering the location of the peak of 
X 4 (t). In order to define a j3 time scale, we first con¬ 
sider an analogous procedure of locating a short-time 
feature in x±{t)- Indeed, as observed in j43|, there 
exists a maximum in X 4 (t) at short times, whose lo¬ 
cation and height depend on the system size, as shown 
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in Fig. m However, we find that the identified time 
scale increases linearly with the length L of the sys¬ 
tem. This behaviour is shown in the inset of Fig[l] 
This plot also shows that this time scale continues to 
increase with increasing L without showing any sign 
of saturation. These observations suggest a phononic 
origin of this time scale. Such a possibility is fur¬ 
ther supported by the fact that this feature disap¬ 
pears when we use Monte Carlo simulations to cal¬ 
culate We opt not to use this procedure for 

identifying Tp. 

In [391 ]. a short time scale Tp was identified as the 
time at which an inflection occurs in a log-log plot of 
the mean squared displacement vs. time. We use the 
same definition to evaluate Tp as a function of system 
size and temperature. In Fig. [2j we show the mean 
squared displacement (MSD) as a function of time in 
a log-log plot and its derivative. The clear dip in the 
derivative and its variation with temperature indicate 
that this estimation of Tp can be done unambiguously 
and without much uncertainty. We use this definition 
of Tp in the present work. 


III. RESULTS 
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FIG. 6 . The relaxation time Tp, plotted as a function 
of the correlation length £, exhibits a power-law depen¬ 
dence, albeit over a small range of values. The line drawn 
through the data points represents a power-law with ex¬ 
ponent z = 0.80. 


In Fig. [3] we show the system-size dependence of 
Tp for different temperatures for the 3dKALJ model. 
For high temperatures the dependence is weak and 
the asymptotic value is reached for small system sizes, 
but this characteristic size becomes increasingly large 
with decreasing temperature. Thus, Tp shows a size 
dependence that is opposite to that displayed by T a 
[341 ] . exhibiting a system-size dependence that may 


normally be expected for a time scale determined by 
an underlying length scale in the system. We note 
that the range over which Tp varies with system size 
and temperature is modest, unlike the behaviour ob¬ 
served for T a . Another important observation about 
the system-size dependence of Tp is that the depen¬ 
dence of this time scale on the system size becomes 
negligibly small if the system size is sufficiently large. 
This is very different from the behaviour of the time 
scale associated with the short-time peak of yqfd) 
shown in Fig. [I] As shown there, this time scale, 
which we attribute to phonon-like excitations, contin¬ 
ues to increase linearly with the length L of the sys¬ 
tem without showing any sign of saturation at large 
values of L. 

We next scale the linear system size L at each tem¬ 
perature by an empirically determined length £(T) 
such that the values of Tp(N,T)/Tp(N —> oo,T) for 
all N and T ( rp{N —> oo, T) is the large -N asymptotic 
value of Tp(N, T)) collapse onto a master curve when 
plotted versus L/£(T). The data collapse obtained 
this way, shown in Fig. [3] is very good, This allows 
for the determination of a /3-regime length scale £(T). 

We then compare the length scale estimated in this 
way with the heterogeneity length estimated in [34 ] 
from quantities obtained at T a . This comparison is 
shown in the inset of the bottom panel of Fig. [3] In 
order to make the comparison, we scale the present es¬ 
timate so that the two length scales match at T = 0.6. 
Length scales obtained by FSS are known only up to a 
multiplicative factor, and hence this procedure does 
not introduce any additional arbitrariness. We find 
that the temperature dependence of the two length 
scales agrees very well. IMCT predicts that £ grows as 
a power law in time up to Tp and stays at the value Tp 
up to T a . Our results therefore confirm this expecta¬ 
tion. Independently of IMCT predictions, the agree¬ 
ment is remarkable in pointing to an intimate connec¬ 
tion between dynamics at short times (/3 regime) and 
long times (a regime). This result suggests that the 
heterogeneity present in the dynamics in the a regime 
(jras been built up already in the /3 regime, and there¬ 
fore in principle, essential information about the a 
relaxation can be obtained by studying short-time dy¬ 
namics. Our result therefore lends support to many 
investigations and ideas that aim to relate dynami¬ 
cal behaviour at short times to long-time structural 
relaxation [24l - l3l[ . 

To see whether the results obtained for the 3dKAL J 
model are generic, we have done similar analysis for 
the 3dIPL and 3dR10 models. In Fig. [2 we have 
shown, in the top left panel, the system-size depen¬ 
dence of the peak value of the four-point susceptibility 
for different temperatures for the 3dIPL model. The 
top right panel shows the results of FSS performed 
for the same data to obtain the dynamic heterogene¬ 
ity length scale. In the bottom left panel we show the 
system-size dependence of the short time scale Tp for 
different temperatures and in the bottom right panel 
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the corresponding FSS of Tp. The scaling collapse 
observed in this case is also very good. We find that 
the temperature dependence of the length scale ob¬ 
tained from the FSS of Tp matches quite well with 
the dynamic heterogeneity length scale, as shown in 
the inset of the bottom right panel of Fig0] Similar 
analysis done for the 3dR10 model confirms that the 
same observations hold for this model too, as shown 
in Fig|5] 

Finally we consider the dependence of rp on the 
extracted length scale £. In Fig. [Gl we show the relax¬ 
ation time as a function of the extracted correlation 
length in a log-log plot. In contrast to T a which ex¬ 
hibits deviations from a power-law dependence on the 
corresponding length scale [H| , we find that a power 
law relation r ~ £ z , holds for Tp. Although such a 
dependence is in qualitative agreement with IMCT, 
we find z ~ 0.80 which is at variance with the IMCT 
prediction^. 


IV. DEPENDENCE ON MICROSCOPIC 
DYNAMICS 

In this section we present results concerning the 
dependence of the time scale Tp on the details of the 
microscopic dynamics. To address this issue system¬ 
atically, we have performed Brownian dynamics (BD) 
simulations using a predictor-corrector scheme [48|. 
This simulation scheme allows one to change system¬ 
atically the friction coefficient to go from the very 
low friction limit (close to a molecular dynamics sim¬ 
ulation) to overdamped dynamics with large friction 
(close to a Monte Carlo simulation). We have changed 
the friction in our simulations over one order of mag¬ 
nitude and studied its effect on the ^-relaxation time 
scale. We find that although a /3 time scale can be 
identified unambiguously for all the values of the fric¬ 
tion coefficients studied, the system-size dependence 
becomes weaker with increasing friction and almost 
completely goes away for the largest value of the fric¬ 
tion coefficient considered here. The system-size de¬ 
pendence of Tp is shown in FigjT] for the 3dKALJ 
model at temperature T = 0.470 for system sizes in 
the range N £ [150,10000]. However, it is important 
to note that the time scale itself remains well-defined 
with increasing friction. This is different from the be¬ 
haviour of the time scale obtained from the short-time 
peak in [Hl : as discussed above, this time scale 
can not be defined for Monte Carlo dynamics because 
the short-time peak in Xi(t) is not present for this 
dynamics. Our results are consistent with previous 
work com par ing the results for different microscopic 
dynamics [49(, which found that the behaviour in the 
early /3 regime is affected by the microscopic dynam¬ 
ics. 

We do not have a full understanding of the ob¬ 
served influence of the microscopic dynamics on the 
system-size dependence of Tp. There are reasons to 



FIG. 7. System-size dependence of the short-time /3- 
relaxation time scale rp at T = 0.470 for the 3dKALJ 
model glass former from Brownian dynamics simulations 
with different values of the friction coefficient. Do indi¬ 
cated in the legend is the microscopic diffusivity of the 
particles and is inversely proportional to the damping co¬ 
efficient. Notice that for small damping (large Do) the 
system-size dependence is very similar to that obtained 
in molecular dynamics simulations. The system-size de¬ 
pendence becomes weaker with increasing damping. The 
time scale for Brownian dynamics is scaled to match the 
long-time mean squared displacement for both molecular 
dynamics and Brownian dynamics. 


believe that the observed behaviour is connected to 
the effects of the properties of the inherent struc¬ 
tures (local minima of the potential energy), whose 
basins of attraction are visited by the system during 
its time evolution, on its dynamics. The dynamics 
in the short-time /3-relaxation regime at low temper¬ 
atures is expected to be strongly influenced by the 
properties of the inherent structures because the sys¬ 
tem should remain confined in the basin of a single 
inherent structure (or in a single metabasin [50], de¬ 
pending on the temperature) during its evolution over 
relatively short times. A recent study [5lJ has shown 
that the dynamics of the 3dKALJ model in the /3 
relaxation regime (up to time scales that are rela¬ 
tively short, but longer than the Tp considered in our 
work), observed in molecular dynamics simulations 
at temperatures near and below the glass transition 
temperature T c of mode coupling theory ( T c ~ 0.435 
for the 3dKALJ model), can be understood from the 
low-energy properties of the relevant inherent struc¬ 
tures. These properties include the eigenvalues and 
eigenvectors of the Hessian matrix evaluated at the 
potential energy minimum (these define the “normal 
modes” of small-amplitude oscillations near the bot¬ 
tom of the basin of an inherent structure) and the 
third and fourth derivatives of the potential energy at 
the minimum (these coefficients determine the effects 
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of anharmonicity on the normal inodes). Thus, the 
system-size dependence of rg obtained in our Newto¬ 
nian molecular dynamics simulations and the length 
scale we have extracted from this dependence should 
be closely related to these properties of the inher¬ 
ent structures. Our observation that the length scale 
extracted from the system-size dependence of rg ob¬ 
tained from molecular dynamics simulations is essen¬ 
tially the same as the length scale of dynamic hetero¬ 
geneity at the a-relaxation time scale then suggests 
that the spatial structure of dynamic heterogeneity 
at time scales of the order of the a-relaxation time is 
closely related to the aforementioned properties of the 
inherent structures. This suggestion receives strong 
support from several experimental [52| and numer¬ 
ical [13, HH, studies that have shown that the 
spatial structure of dynamic heterogeneity at the a- 
relaxation time scale is closely related to the structure 
of the eigenvectors associated with some of the low- 
lying eigenvalues of the Hessian matrix evaluated at 
the appropriate inherent structures. All these obser¬ 
vations suggest that the short-time (/3-relaxation) and 
long-time (a-relaxation) dynamics observed in molec¬ 
ular dynamics simulations are closely related to each 
other through the low-energy properties of the rele¬ 
vant inherent structures. In Brownian dynamics sim¬ 
ulations with a large friction coefficient, the influence 
of the properties of the inherent structures on the 
short-time dynamics is masked by the strong damping 
of the normal modes. This may reduce or completely 
eliminate some of the effects of the structural proper¬ 
ties of the inherent structures on the dynamics. The 
system-size dependence of rg which, as argued above, 
is closely related to the low-energy properties of the 
inherent structures, may be one of the effects that 
disappear in the presence of strong friction. These 
arguments provide a rationalization of the observed 
effects of strong friction on the system-size depen¬ 
dence of rg obtained from Brownian dynamics simu¬ 
lations and suggest that the low-energy properties of 
the inherent structures visited by the system during 
its time evolution play an important role in both /3 
and a relaxation processes. 


V. DISCUSSION 


In this work, we have studied the system-size de¬ 
pendence of the /3 relaxation time for three model 


liquids in three dimensions. In each case, we find 
that the time scale initially increases with increasing 
system size and saturates for large values of the sys¬ 
tem size, exhibiting behaviour conforming to usual 
expectations for the size dependence of a quantity 
that depends on a length scale. This is unlike the 
system-size dependence of the a relaxation time that 
is found to decrease with increasing system size 341. 
On the other hand, the length scale extracted from 
FSS of this /3 relaxation time scale matches very well 
with the length scale extracted from the FSS of the 
four-point susceptibility Xa (T) measured at the a- 
relaxation time scale. This is in agreement with ex¬ 
pectations based on inhomogeneous mode coupling 
theory, as also our observation that the asymptotic 
value of the /3 time scale, obtained for large systems, 
exhibits a power-law dependence on the length scale. 
The value of exponent of this dependence, however, 
does not match IMCT predictions (^ ~ 5.4) )3j. 
Our results are also consistent with ideas that re¬ 
late dynamics at short times to long-time structural 
relaxation [234^ . We have also presented intriguing 
results for the dependence of the short-time dynamics 
on the microscopic dynamics used in the simulation 
and suggested a possible explanation of the observed 
behaviour in terms of the role played by certain low- 
energy properties of the relevant inherent structures 
(i.e. inherent structures whose basins are visited by 
the system during its time evolution) in the dynamics 
of the system. 


The results reported here have a strong connection 
with those in [55j |. in which oscillatory shear simu¬ 
lations were performed to study the effects of short- 
time /3-relaxation in a supercooled two-dimensional 
glass-forming liquid at both high and low tempera¬ 
tures. The loss modulus measured in this study was 
shown to be related to /3-relaxation and the effect of 
this relaxation was found to decrease sharply by the 
introduction of a small amount of pinning disorder in 
the system. The cooperative displacements of indi¬ 
vidual particles over the /3-relaxation time scale were 
found to diminish very rapidly with the introduction 
of pinning disorder. Based on these observations, it 
was concluded that /3-relaxation is also cooperative 
in nature, similar to a-relaxation. Our study clearly 
shows that ^-relaxation is indeed cooperative and the 
cooperativity is likely to have the same origin as that 
in a-relaxation. 
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